Dynamics of Energy Fluctuations in Equilibrating and Driven-Dissipative Systems 
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When two isolated system are brought in contact, they relax to equihbrium via energy exchange. 
In another setting, when one of the systems is driven and the other is large, the first system reaches 
a steady-state which is not described by the Gibbs distribution. Here, we derive expressions for 
the size of energy fluctuations as a function of time in both settings, assuming that the process is 
composed of many small steps of energy exchange. In both cases the results depend only on the 
average energy flows in the system, independent of any other microscopic detail. In the steady-state 
we also derive an expression relating three key properties: the relaxation time of the system, the 
energy injection rate, and the size of the fluctuations. 
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FIG. 1. The energy fluctuations in the two set-ups. (a) Ex- 
change of energy between two systems, (b) A system driven 
by an external force and attached to a bath, (c) A typical 
evolution of E\ (t) (dashed line) fluctuates around the aver- 
age (El) (t) (solid line). Eq. (O relates the size of these 
fluctuations to the average (Ei) (t). 



In this paper we consider two closely related non- 
equilibrium problems. In the first problem, two systems 
which are coupled to each other but isolated otherwise, 
are allowed to exchange energy, see Fig. [T{a). The sys- 
tems start with arbitrary initial energies and eventually 
reach equilibrium. It is natural to ask: How do the ini- 
tial energies evolve in time as the two systems approach 
equilibrium? For example, one might imagine measuring 
the energy of a tea cup as it cools, or the equilibration 
of a mesoscopic system of two atomic gases, initially pre- 
pared at two different temperatures. In the second prob- 
lem, one of the two systems is also driven by an external 
protocol, see Fig. [IJb). This is achieved, for example, by 
applying a time-varying field which repeatedly returns 
to its initial form. When the second system is much 
larger than the first, it acts as a dissipative bath, and 
the first system eventually settles to a non-equilibrium 
steady-state. This scenario serves as a generic model for 
driven-dissipative systems, which describe a broad range 
of phenomena [H-Q. Here, one can ask how the first sys- 
tem reaches a steady-state, and what are the properties 
of this non-Gibbsian steady-state. 



As the dynamics of a system are affected by the de- 
tailed microscopic state, repeating the same experiment 
will lead to different outcomes. Specifically, a measure- 
ment of the energy as a function of time will yield differ- 
ent results, see Fig. [Ijc). The variations between exper- 
iments might average-out in large, thermodynamic sys- 
tems, or when the driving protocol applied is quasi-static. 
However, they can be significant when the drive is not 
quasi-static and in small or mesoscopic systems, which 
are of current experimental interest [^-l6| . Here we quan- 
tify these energy fluctuations by studying the variance of 
the energy measurements in repeated experiments. The 
dependence of these fluctuations on the dynamics makes 
general statements scarce, and one typically has to resort 
to the study of speciflc models. 

In this Letter we show that when the changes in energy 
are small and slow (but still irreversible), general state- 
ments about the energy fluctuations can be made. The 
results are insensitive to almost all microscopic details of 
the systems, depending only on the average energy flows 
from the drive to the system and between the systems 
as a function of time, and on the density of states. We 
stress that the assumptions made do not imply that the 
combined system (composed of systems 1 and 2) is close 
to equilibrium, but only that each of the systems sepa- 
rately is close to equilibrium within its energy shell. Our 
main results are: (1) Eqs. ^ and ([7]), which quantify 
the variance of the energy fluctuations as the system ap- 
proaches its steady-state (which is equilibrium when no 
drive is present). (2) Eq. which relates three main 
quantities at the steady-state: the variance of the energy 
fluctuations, the average rate of energy flow through the 
system, and the relaxation time of energy fluctuations. 
The validity of the results is illustrated in a system of 
colliding hard spheres. 

To derive them we consider the evolution of the ener- 
gies in the {Ei , E2 ) plane, where Ei , E2 are the energies 
of systems 1 and 2 respectively. Consider a series of small 
changes in the energies, each taking place over a time 
interval At. We assume that Tfj^ ^ At, where is the 
relaxation time of each of the isolated systems separately. 
When this time scale separation holds, the statistics of 
the energy changes AEi , AE2 during the time interval 
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from t to t + At depend only on the energies {Ei,E2) 
at time t. The time evolution of the probability distri- 
bution Pi2 {El, E2) is then governed by a Fokker-Planck 
equation 
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dtPl2 ^-Y.^E^ (^^^12) + \Y. (B.,-Pi2) (1) 

i=l i-/j=l 

where Ai, A2, i?i2, B21, ^22 are all functions of 
{Ei,E2), and B12 ~ -621- These function are related 
to the first two moments of the changes in Ei , E2 during 
a short time At 0,11]: 
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The equation is valid when higher cumulants, e.g. 
(AEiAEjAEk) /At, are small compared to the Ai and 
Bij functions. 

In both scenarios, of equilibrating systems and driven- 
dissipative systems, one can take the Ai and Bij func- 
tions to depend on Ei only. For equilibrating systems, 
this is possible when the initial total energy Etotai — 
El + E2 is fixed, so that E2 can be considered to be a 
function of Ei. In the case of driven-dissipative systems, 
E2 drops completely from the equations when we take 
system 2 to be much larger than system 1 . As shown be- 
low, this is because system 2 acts as a thermal bath whose 
properties are insensitive to the changes in E2. It is then 
more convenient to work with the marginal probability 
distribution of Ei alone: P {Ei) = f dE2Pi2 {Ei, E2). 
Integrating Eq. ([TJ over E2 we find 
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dtP = -Oe, (AiP) + -di,^ {BiiP) 



(2) 



using P12 {E2 — ±00) = 0. While only the functions 
Ai and Bn appear in this equation, the interaction with 
system 2 still affects the energy of system 1, via the forms 
of the functions Ai, Bn. This is contained in the relation 
which is derived below 



2Ai^ 2/32/ PiAp^iPi- 132) Bii 



(3) 



where Ap = Ai + A2 = dt {Etotai) is the rate of en- 
ergy injected into the system by the drive. The inverse 
temperatures are defined by /3i {Ei) — dsiSi {Ei), and 
/32 {E2) = dE2S2 {E2), where 5i^2 are the (microcanoni- 
cal) entropies of systems 1,2 respectively. /3i and /32 are 
well-defined functions, depending only on the density of 
states of the system, and unrelated to the driving mecha- 
nism and the interaction between the systems. Moreover, 
/3i can be very different from /32. Eq. (jSj is ultimately 
based on Liouville's equation, or the unitarity of the dy- 
namics in quantum cases. In the driven case we also 
assume that the energy flow from the drive and between 
the systems are statistically independent processes, see 
discussion below. Eq. is exact up to corrections of 
order 1/N , where N is the number of degrees of freedom 



of the smaller of the two systems. In the case of equili- 
brating systems = 0, and the relation Eq. ^ reduces 
to 2Ai = (/3i — (32) Bii. Here, as expected, on average 
energy flows from high to low temperatures. 

The drive is implemented by varying the Hamiltonian 
of system 1 in time (e.g., by applying a time- varying ex- 
ternal field). We consider drives where the Hamiltonian 
repeatedly returns to its original form (i.e., an oscillat- 
ing field). At the steady-state, when the Hamiltonian 
is changed adiabatically, returning to the original form 
leaves the energy of the combined system unchanged. 
Thus, the changes in the energy will only be due to irre- 
versible effects. 

Before deriving Eq. ([3]) we consider several of its conse- 
quences in the two scenarios, of equilibrating and driven- 
dissipative systems. Wherever possible, we present the 
results in a unified way where the case of equilibrating 
systems is obtained by setting Ap = 0. 

Approach to steady-state - We start by considering the 
approach of the combined system (composed of systems 
1 and 2) to its steady-state. If no driving is present (sce- 
nario 1), this steady-state is thermal equilibrium. We 
derive an expression for the evolution of the variance 
af = — {Ei)'^ during the entire equilibration pro- 

cess. Proceeding similarly to @, we take the first two 
moments with respect to Ei of Eq. ^ 



dt {El) = {Ai) , 

dt(Tl = {Bii)+2{{AiEi) 



{M) {El)) 



(4) 



If the distribution is narrow enough (valid up to 1/iV 
corrections, see discussion after Eq. ([7])), {Ai) can be 
assumed to depend on {Ei) alone, and the change in {Ei) 
will be monotonic. Combining the two equalities in Eq. 
(HI) and linearizing Ai within the width of the probability 
distribution, we find 



d{Ei 



2Z{{Ei)) + 2- 



{Ai 



^tHEi)) , (5) 



where Z {{Ei)) = Buj (2Ai). Solving the ordinary dif- 
ferential equation Eq. ([5]) and using Eq. ([3]) we find for 
the equilibrating systems that the variance is given by 

2 2 A\{{Ei)) 

<^l.eq {{El)) =(Ji„- 
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2Al{{Ei)) 



(El) 



^E,)„ AI {E') {fil - h) 
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Here (i?i)o and are (£'1) and cr\ respectively at the 
initial time. Recall that Etotai is held constant in this 
expression. It is easy to extend these results when Etotai 
varies between experiments. It is interesting to note that 
this expression is identical to that obtained for a single 
driven isolated system Q when /32 is set to zero. This 
means that within this theory, driving a system is for- 
mally equivalent to attaching it to a bath with infinite 
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temperature. It is straightforward to show, that when 
system 2 is a bath, so that ^2 can be taken to be a con- 
stant, the width crj approaches the equihbrium value: 



ksT'^C, were E^„ is the equilibrium value 



of {El), and C is the heat-capacity (see e.g., [lOj). To 
see this, note that at equilibrium Ai must vanish, and 
/3i = h- Therefore the entire expression for a\ is con- 
where 



trolled by the final approach of E to Eeq 



Pl~P2 + {dE,P2)E^, (E^ 

{E-E,q) 



Note that away 
[E) need not be 



and the equilibrium expression follows, 
from the final equilibration regime Ai 
linear. 

In the case of driven-dissipative systems (when system 
2 is large), we obtain for the variance 
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where Z = [I - ^^Ap I {^iAi)\ / {f3i - (32). Eqs. © and 
([7|) are our main results for the approach to steady-state. 
They predict the size of fluctuations in Ei around its av- 
erage value. They depend only on the rates of energy 
injection into the system Ap (which is zero for equili- 
brating systems) and the rate of energy transfer to the 
bath A2. In principle both these quantities can be mea- 
sured separately. Ap can be measured by the rate of 
energy absorption when system 1 is isolated, and A2 in 
an equilibration experiment without the drive. This is 
a consequence of our assumption of statistical indepen- 
dence of the driving and the mechanism of interaction 
between the systems. In addition we comment that Eqs. 

and d?]) imply that crf ^^^ and cr^ scale as E oc N 
when A is a homogeneous function of N (e.g., extensive 
in N). This justifies self-consistently our assumption on 
the narrowness of the distribution. 

Steady-state fluctuations - The framework described 
above can also be used to study fluctuations in the 
steady-state of driven-dissipative systems, specifically 
fluctuations of Ei around {Ei). At the steady-state the 
probability distribution Pg {Ei) is independent of time. 
Using dt {El) = {A) = A{{Ei)), and noting that at the 
steady-state A{{Ei)) must vanish, we expand Ai and 
Bii to lowest order in ei = Ei — Ei 
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Al = ei 

T 



Bii — Bs 



(8) 



where Bs and r are constants. Equivalently, in this 
regime the Fokker-Planck equation describes the Brown- 
ian motion of the energy in a harmonic potential ei = 
— ci/t -I- -v/Bjry, where the white noise 77 (i) satisfles 
(r) {t) rj {t')) ^ 5{t — t'). T is then interpreted as the re- 
laxation time, as can be seen from the two time correla- 
tion function 

(ei {ti)ei {t2))^^i 



H«2-tl|/r 



The variance of the energy fluctuations is given by af = 
ei {ti)^) ^ Bst/2. 

When Ap = 0, Eqs. © and © imply that d = 
— ^-§^ (/3i — ^2)- Then expanding /3i around as done 

above we flnd that = Bst/2 = — (j)E^Pij which 

again reproduces the canonical distribution width. The 
present derivation gives a dynamic interpretation to this 
formula. 

When Ap ^ namely for a driven-dissipative system 
we flnd, using Ai (£'?) = in Eq. ©, and af = Bst/2, 
that 



tAi 



P2 



(/32-/3i)a? 



(9) 



This is our main result for the steady-state of driven- 
dissipative systems. Ap is the rate of energy injected 
to the system from the drive. In the steady-state, this 
energy is then dissipated into the bath. This expression 
therefore relates three central quantities characterizing 
the steady-state: the size of the energy fluctuations cr^, 
the rate of energy dissipation Ap, and r which is the 
relaxation time in the steady-state [T]| . 

MD Simulations - Before proving the key relation Eq. 

, we illustrate our main results on a gas of hard-sphere 
particles in a box, simulated by an event-driven molec- 
ular dynamics simulation p^ . The gas is composed of 
A^i particles of mass mi and N2 particles of mass TO2 , all 
of equal size, corresponding to systems 1 and 2 respec- 
tively. Although the entropy of the two systems between 
collisions indeed factorizes, the collision process involves 
a strong interaction, which changes the velocities of the 
particles by a significant amount. A collision calcula- 
tion shows that if the two masses are very different, the 
energy transfer in each collision is small. In this case 
energy transfer occurs over many collisions, fulfilling the 
assumption of time-scale separation (see above) . In what 
follows we take mi = 10~^, 7712 = 1. (Throughout we use 
arbitrary units). The box is a unit cube with reflecting 
boundary-conditions, and the particles are taken to oc- 
cupy a volume fraction of 0.05. 

We first consider the approach to equilibrium of two 
systems in contact. Fig. [TJa), to be compared with the 
predictions of Eq. ©. We take Ni = 30 for the first 
systems and A^2 — 20 for the second system. A''i , A^2 are 
chosen to be relatively small in order to test the theory on 
a mesoscopic system. The initial velocities are sampled 
from a Maxwell-Boltzmann distribution with /3i = 60 
and 132 = 3, corresponding to average energies per parti- 
cle of {El) /Ni = 0.025 and {E2) /N2 = 0.5. We start aU 
runs from a fixed total energy Etotai — {Ei) + {E2), by 
preforming a (small) rescaling of the m2-particles' veloc- 
ities. Gathering statistics over many runs, we calculate 
at each time the average energy {Ei) {t) and the vari- 
ance a\ {t). The function Ai {{Ei)) is obtained by plot- 
ting Al {t) = d {El) /dt as a function of {Ei) {t). Given 
Ai{{Ei)) [Hi we use Eq. 1^ to predict af{{Ei)), and 
find a good fit with the simulation results, see Fig. [51 
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FIG. 2. Equilibration of a system of 50 particles with two 
different masses. Plotted are simulation results for vs. 
(dots) compared to the theoretical prediction (solid line) , Eq. 
Q. Inset: Example of -Bi (t) in a single run (solid line), and 
the average energy (t) (dashed line), used to calculate 
Ai. 
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FIG. 3. Test of Eq. ^ for a driven-dissipative system. 



We now turn to the driven-dissipative scenario, and 
test Eq. dni) for the steady-state. We run simulations on 
a system with A^i = 10 and N2 = 50 particles. System 1 
is driven by applying short impulses to the mi-particles, 
changing their velocity by Av — FSt/mi, where F is a 
constant force and St the impulse duration. This impulse 
is applied at a constant rate. In order to mimic the be- 
havior of a very large system 2, the velocity of the TO2 
particles is changed upon reflection from the wall [14] , so 
as to maintain a constant (£'2). The quantities t,Af,/3i 
and a I are computed from the numerics. The energy of 
system 2 is maintained at {E2) /N2 — 1/2, or (3m = 3. 
Fig. [3] shows the results obtained for the two sides of 
Eq. as a function of Pi/ ^2 for different strengths of 
the drive. Good agreement is found over a wide range 
of drive strengths, and temperature differences. In this 
range, Ap increases by a factor of 1000, al by a factor of 
17, and the relaxation time r decreases by a factor of 3. 

Derivation of the eq. ^ - We end by deriving Eq. ^ . 
To do so we look at two times t, t+At in which the driving 



protocol has returned to its original state, i.e. where 
H (t) = H {t + At) , where H is the Hamiltonian of the 
combined system. As stated before, we assume that both 
subsystems are relaxed in their respective energy shells, 
with energies Ei and £2- In particular, this requires that 
At Tfl. We denote the changes in Ei,E2 during the 
time interval At by AEi , AE2 respectively, and define 
AEb and AEp via 

AE'i = AEp - AEb , AE2 = AEb ■ 

AEb is the energy transferred from system 1 to system 
2, and AEp is the work done on system 1 by the external 
drive. 

Under these assumptions, Liouville's theorem or 
the unitarity of the dynamics, together with micro- 
reversibility of the dynamics, imply a Crooks relation for 
the combined isolated system [15l-[l7| 



= Pe,+ae,,e.+ae. i-AEi, -AE2) eS.(^i+A^^.)+s.(i^.+A^..) 

where Pei,E2 {AEi, AE2) is the probability of a transi- 
tion {Ei,E2) (El + AEi,E2 + AE2), and P is defined 
similarly, only with respect to the the reversed protocol, 
defined by the dynamics generated by the time-reversed 
Hamihonian, H {t') ^ H {t + At - t'). 

We assume that: (1) AEi,AE2 are small so that the 
transition probabilities depend weakly on the initial en- 
ergies Pei+aEi,E2±AE2 = Pei,E2j with corrections which 
are of order Q- (2) AEb and AEp are statistically 
independent quantities. This happens when the interac- 
tion with the bath is independent from the driving pro- 
cess, e.g., when the drive and interaction processes act 
on different modes, on different parts of the system, at 
different times, etc.. These imply that 

{AEp, 0) P^ {-AEb, AEb) e-'^'^^^-'^'^'-^'^^^^ 

= P^ {AEp,0)P^ {-AEb,AEb) . 

Integrating over AEp, AEb gives a Jarzynski re- 



lation (e 



3iAEf 



-{P2-Pi)AE, 



) = 1. Rear- 
again integrating gives 



PiAEf 



) 



ranging the equation and 

^g-ftAB^^ ^ (^^-(P2-h)AEs^ yielding (e" 

^g-(02-/3i)A£;i3 ^ = 1. The second relation is the exchange 
fluctuation relation [l^ [3| ; The first is a variant of the 
Jarzynski relation for isolated systems Expanding 
both to second order we find 



2 (AEp) = /3i (AE],) , 2 (AEb) = {P2 - Pi) {AEI) . 

From independence it follows that {AEpAEi) = (Ai?|,), 
or (A£;|) = (A^;^) - (AEj,). The second equation be- 
comes 2 (AEb) = {P2 - Pi) ({^EI) - (AEj,)), so that 

2 (AEi) - 2P2/P1 {AEp) = {Pi - P2) {AEl) . 

After dividing by At and using the definitions of the 
Ai,Bi quantities, yields Eq. ([3]). 
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